Usenet 1993 July/InfoMagic USENET CD-ROM July 1993.ISO/sources/unix/volume21/rayshade/part03

Internet Message Format | 1990-03-21 | 56.1 KB

Subject:  v21i010:  A ray tracing program, Part03/08
Newsgroups: comp.sources.unix
Sender: sources
Approved: [e-mail]

Submitted-by: Craig Kolb <[e-mail]>
Posting-number: Volume 21, Issue 10
Archive-name: rayshade/part03

#! /bin/sh
# This is a shell archive.  Remove anything before this line, then unpack
# it by saving it into a file and typing "sh file".  To overwrite existing
# files, type "sh file -c".  You can also feed this as standard input via
# unshar, or by typing "sh <file", e.g..  If this archive is complete, you
# will see the following message at the end:
#        "End of archive 3 (of 8)."
# Contents:  BLURB.UTAH src/bounds.c src/cone.c src/cylinder.c
#   src/intersect.c src/light.c src/object.c src/ray_options.c
#   src/voxels.c
PATH=/bin:/usr/bin:/usr/ucb ; export PATH
if test -f 'BLURB.UTAH' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'BLURB.UTAH'\"
else
echo shar: Extracting \"'BLURB.UTAH'\" \(6402 characters\)
sed "s/^X//" >'BLURB.UTAH' <<'END_OF_FILE'
X
X               THE UTAH RASTER TOOLKIT
X
XThe Utah Raster toolkit is a collection of programs and C routines for
Xdealing with raster images commonly encountered in computer graphics.  It
Xprovides the following major functions:
X
X    * A device and system independent image format for storing images
X      and information about them.  Called the RLE format, it uses
X      run length encoding to reduce storage space for most images.
X
X    * A library of C routines for reading, writing and manipulating
X      images stored in the RLE format.
X
X    * A collections of programs for manipulating and displaying RLE
X      images.
X
X
XThe Format:
X
X  The device independent RLE file has two parts, the header, which stores
X  information about the image (size, position, channel information,
X  color maps, comments, etc), and the actual image data in a run length
X  encoded format.  The RLE format often requires about a third of the
X  available space for most "image synthesis" style images.  If the image
X  does not compress well, the RLE format stores it as straight pixel data
X  with little extra overhead.  The format has been developed over the past
X  five years at Utah.
X
XThe Library:
X
X  C routines are provided for setting up and reading the image header,
X  and for reading and writing the image a scanline at a time.  Images can
X  be read or written using two different methods.  Using the "row" method,
X  the library performs the RLE encoding and decoding.  With the "raw" method,
X  scanlines are constructed directly with RLE opcodes.  Additional routines
X  are available for generating dither matrices (e.g., for display programs
X  running on devices with less than 24 bits of color).
X
XThe Tools:
X  applymap   - Apply color map values to pixel values.
X  avg4          - Downfilter an image by 1/4, generating a matte channel if one
X           didn't previously exist
X  comp          - Digital image compositor.  Provides the operations over, atop,
X           in, out, xor, plus, minus and diff on two images.
X  crop          - Crop an image.
X  dvi2rle    - Convert TeX output into anti-aliased images.
X  fant          - Rotate and/or scale in image by an arbitrary (float) value.
X  mcut       - Quantize an image from 24 to eight bits using the median cut
X           algorithm.
X  mergechan  - Merge several channels from different files into a single
X           RLE file.
X  pyrmask    - Blend images using Gaussian pyrimids.
X  repos      - Change the position in the RLE header.
X  rleClock   - Generate an image of a clock.
X  rleaddcom  - Add a comment to the RLE file's header.
X  rlebg      - Generate a solid or variable background.
X  rlebox     - Find the actual non-background area of an image.
X  rleflip    - Rotate an image by 90/180 degree increments.
X  rlehdr     - Dump the contents of the RLE header in human readable form.
X  rlehisto   - Generate the histogram of an RLE file.
X  rleldmap   - Load a color map into an RLE file from a variety of sources.
X  rlemandl   - Generate Mandlebrot sets as RLE files.
X  rlenoise   - Adds random noise to an image.
X  rlepatch   - Overlay several smaller images over a larger one.
X  rlescale   - Generates gray scale and color scale RLE files.
X  rlesetbg   - Set the background color stored in the RLE header.
X  rlesplit   - Split a file containing several images into several files.
X  rleswap    - Swap, copy or delete channels in an RLE file.
X  rletops    - Convert an RLE image to PostScript (graylevel).
X  rlezoom    - Enlarge an image with pixel replication.
X  smush      - Perform a simple Gaussian filter on an image.
X  to8          - Convert a 24 bit RGB image to an eight bit dithered one.
X  tobw          - Convert 24 bits to 8 bits black and white.
X  unexp      - Convert an "exponential" image to a displayable one.
X  unslice    - Quickly assemble an image from several horizontal strips
X
X  Format conversion programs are provided for:
X       - Simple pixel streams (color & B&W)
X    - Targa image format
X    - Cubicomp image format
X    - PostScript
X    - MacPaint
X    - Sun rasterfiles
X    - Wastatch paint systems
X
X  Display programs are provided for:
X  getap         - Apollo workstations
X  getbob     - HP Series 300 ("bobocat") running Windows 9000
X  getcx3d    - Chromatics CX1500 display
X  getfb      - BRL "libfb" displays
X  getgmr     - Grinnell GMR-27 (remember those?)
X  getX         - Workstations running the X window system
X  getX11     - Workstations running X11
X  getOrion   - Orion displays
X  getren     - HP 98721 "Rennasance" display
X  getsun     - Suns running Suntools
X  getmac     - Macintosh.
X  getmex     - Iris running Mex
X  getqcr     - Photograph images with the Matrix QCR-Z camera.
X  getiris    - Iris in raw 24 bit mode.
X             - [Note display programs for a particular device are
X        simple to add]
X
X  All the tools are designed to pipe together, so they can be used as
X  filters on images much like the standard Unix tools filter text.
X
XPlus:
X
X  The raster toolkit also includes Unix man pages for the library and
X  commands, some sample images, and additional documentation.
X
XSystem Requirements:
X
X  We have successfully ported the Raster Toolkit to a number of Unix
X  systems, including 4.2/4.3bsd (Vax, Sun, etc), Apollo Domain/IX, HP
X  Series 3000, SGI Iris, Gould UTX.  Display programs are included for
X  several devices.  Creating display programs for additional devices is
X  a straightforward task.
X
XDistribution:
X
X  For ARPAnet sites, the toolkit may be obtained via anonymous FTP to the
X  site cs.utah.edu, in the file pub/toolkit-2.0.tar (or, if you cannot FTP
X  that large a file at once, in pub/toolkit-2.0.tar.1, pub/toolkit-2.0.tar.2
X  and pub/toolkit-2.0.tar.3).  Sites not on the ARPAnet can obtain the Raster
X  Toolkit on a 9-track, 1600 bpi tar format tape by sending check or
X  money order for $200.00, payable to the Department of Computer Science,
X  to:
X
X    Attn: Utah Raster Toolkit, Loretta Looser
X    Department of Computer Science
X    University of Utah
X     Salt Lake City, UT, 84112
X
X  Courtesy Mike Muuss at BRL, the Raster Toolkit is also included as
X  contributed software in the BRL-CAD distribution.
X
X  [Note: because of the size of the distribution, we can not distribute
X  it via mail or UUCP]
X
X  Although the Raster Toolkit software is copyrighted, it may be freely
X  re-distributed on a "GNU-like" basis.
X
XFor further technical information on the Raster Toolkit, send mail
Xto:
X    [e-mail]        (ARPA)
X    {ihnp4,decvax}!utah-cs!toolkit-request    (UUCP)
X
X
END_OF_FILE
if test 6402 -ne `wc -c <'BLURB.UTAH'`; then
    echo shar: \"'BLURB.UTAH'\" unpacked with wrong size!
fi
chmod +x 'BLURB.UTAH'
# end of 'BLURB.UTAH'
fi
if test -f 'src/bounds.c' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'src/bounds.c'\"
else
echo shar: Extracting \"'src/bounds.c'\" \(5679 characters\)
sed "s/^X//" >'src/bounds.c' <<'END_OF_FILE'
X/*
X * bounds.c
X *
X * Copyright (C) 1989, Craig E. Kolb
X *
X * This software may be freely copied, modified, and redistributed,
X * provided that this copyright notice is preserved on all copies.
X *
X * There is no warranty or other guarantee of fitness for this software,
X * it is provided solely .  Bug reports or fixes may be sent
X * to the author, who may or may not act on them as he desires.
X *
X * You may not include this software in a program or other software product
X * without supplying the source, or without informing the end-user that the
X * source is available for no extra charge.
X *
X * If you modify this software, you should include a notice giving the
X * name of the person performing the modification, the date of modification,
X * and the reason for such modification.
X *
X * $Id: bounds.c,v 3.0 89/10/27 02:05:46 craig Exp $
X *
X * $Log:    bounds.c,v $
X * Revision 3.0  89/10/27  02:05:46  craig
X * Baseline for first official release.
X *
X */
X#include <stdio.h>
X#include <math.h>
X#include "constants.h"
X#include "typedefs.h"
X#include "funcdefs.h"
X
X#ifndef HUGE
X#define HUGE    1.701411e38
X#endif
X
X/*
X * Ray-bounding box intersection test.
X */
Xdouble
XIntBounds(ray, bounds)
XRay *ray;
Xdouble bounds[2][3];
X{
X    double t, tmin, tmax, bmin, bmax;
X    double dir, pos;
X    extern unsigned long BVTests;
X
X    BVTests++;
X    tmin = 0.;
X    tmax = FAR_AWAY;
X
X    dir = ray->dir.x;
X    pos = ray->pos.x;
X
X    if (dir < 0) {
X        bmin = bounds[HIGH][X];
X        bmax = bounds[LOW][X];
X    } else {
X        bmax = bounds[HIGH][X];
X        bmin = bounds[LOW][X];
X    }
X
X    if (dir != 0.) {        /* check x-faces */
X        t = (bmax - pos) / dir;
X        if (t < 0.)
X            return 0.;
X        if (t <= tmax)
X            tmax = t;
X        t = (bmin - pos) / dir;
X        if (t >= 0.) {
X            if (t > tmax)
X                return 0.;
X            tmin = t;
X        }
X    } else if (pos < bmin || pos > bmax)
X        return 0.;
X
X    dir = ray->dir.y;
X    pos = ray->pos.y;
X
X    if (dir < 0) {
X        bmin = bounds[HIGH][Y];
X        bmax = bounds[LOW][Y];
X    } else {
X        bmax = bounds[HIGH][Y];
X        bmin = bounds[LOW][Y];
X    }
X
X    if (dir != 0.) {        /* check y-faces */
X        t = (bmax - pos) / dir;
X        if (t < 0.)
X            return 0.;
X        if (t <= tmax) {
X            if (t < tmin)
X                return 0.;
X            tmax = t;
X        }
X        t = (bmin - pos) / dir;
X        if (t >= tmin) {
X            if (t > tmax)
X                return 0.;
X            tmin = t;
X        }
X    } else if (pos < bmin || pos > bmax)
X        return 0.;
X
X    dir = ray->dir.z;
X    pos = ray->pos.z;
X
X    if (dir < 0) {
X        bmin = bounds[HIGH][Z];
X        bmax = bounds[LOW][Z];
X    } else {
X        bmax = bounds[HIGH][Z];
X        bmin = bounds[LOW][Z];
X    }
X
X    if (dir != 0.) {        /* check z-faces */
X        t = (bmax - pos) / dir;
X        if (t < 0.)
X            return 0.;
X        if (t <= tmax) {
X            if (t < tmin)
X                return 0.;
X            tmax = t;
X        }
X        t = (bmin - pos) / dir;
X        if (t >= tmin) {
X            if (t > tmax)
X                return 0.;
X            tmin = t;
X        }
X    } else if (pos < bmin || pos > bmax)
X        return 0.;
X
X    return tmin;
X}
X
X/*
X * Transform an object's bounding box by the given transformation
X * matrix.
X */
Xtransform_bounds(trans, objbounds)
XTransInfo *trans;
Xdouble objbounds[2][3];
X{
X    Vector v, tmp;
X    double bounds[2][3];
X    int x, y, z;
X
X    init_bounds(bounds);
X
X    /*
X     * Find bounding box of transformed corners of bounding box.
X     */
X    for (x = 0 ; x < 2; x++) {
X        v.x = objbounds[x][X];
X        for (y = 0; y < 2; y++) {
X            v.y = objbounds[y][Y];
X            for (z = 0; z < 2; z++) {
X                v.z = objbounds[z][Z];
X                tmp = v;
X                transform_point(&tmp, trans);
X                if (tmp.x < bounds[LOW][X])
X                    bounds[LOW][X] = tmp.x;
X                if (tmp.x > bounds[HIGH][X])
X                    bounds[HIGH][X] = tmp.x;
X                if (tmp.y < bounds[LOW][Y])
X                    bounds[LOW][Y] = tmp.y;
X                if (tmp.y > bounds[HIGH][Y])
X                    bounds[HIGH][Y] = tmp.y;
X                if (tmp.z < bounds[LOW][Z])
X                    bounds[LOW][Z] = tmp.z;
X                if (tmp.z > bounds[HIGH][Z])
X                    bounds[HIGH][Z] = tmp.z;
X            }
X        }
X    }
X
X    for (x = 0; x < 3; x++) {
X        objbounds[LOW][x] = bounds[LOW][x];
X        objbounds[HIGH][x] = bounds[HIGH][x];
X    }
X}
X
Xinit_bounds(bounds)
Xdouble bounds[2][3];
X{
X    bounds[LOW][X] = bounds[LOW][Y] = bounds[LOW][Z] = HUGE;
X    bounds[HIGH][X] = bounds[HIGH][Y] = bounds[HIGH][Z] = -HUGE;
X}
X
X/*
X * Walk through a linked-list of objects.  If the object is unbounded,
X * unlink it it from the list and add it to the 'unbounded' list.
X * If the object is bounded, enlarge the given bounding box if
X * necessary.  Return pointer to unbounded list.
X */
XObjList *
Xfind_bounds(list, bounds)
XObjList **list;
Xdouble bounds[2][3];
X{
X    ObjList *ltmp, *prev, *oltmp;
X    ObjList *unbounded;
X    Object *otmp;
X
X    init_bounds(bounds);
X    prev = unbounded = (ObjList *)0;
X
X    for (ltmp = *list; ltmp; ltmp = ltmp->next) {
X        otmp = ltmp->data;
X        if (otmp->bounds[LOW][X] > otmp->bounds[HIGH][X]) {
X            /*
X             * Object is unbounded -- unlink it...
X             */
X            if (prev)
X                prev->next = ltmp->next;
X            else
X                *list = ltmp->next;
X            /*
X             * And add it to unbounded object list.
X             */
X            oltmp = (ObjList *)Malloc(sizeof(ObjList));
X            oltmp->data = otmp;
X            oltmp->next = unbounded;
X            unbounded = oltmp;
X        } else {
X            /*
X             * Object is bounded.
X             */
X            enlarge_bounds(bounds, otmp->bounds);
X            prev = ltmp;
X        }
X    }
X    return unbounded;
X}
X
X#define SetIfLess(a, b)        (a = (a) < (b) ? (a) : (b))
X#define SetIfGreater(a, b)    (a = (a) > (b) ? (a) : (b))
X
X/*
X * Find bounding box of the union of two bounding boxes.
X */
Xenlarge_bounds(old, new)
Xdouble old[2][3], new[2][3];
X{
X    SetIfLess(old[LOW][X], new[LOW][X]);
X    SetIfLess(old[LOW][Y], new[LOW][Y]);
X    SetIfLess(old[LOW][Z], new[LOW][Z]);
X    SetIfGreater(old[HIGH][X], new[HIGH][X]);
X    SetIfGreater(old[HIGH][Y], new[HIGH][Y]);
X    SetIfGreater(old[HIGH][Z], new[HIGH][Z]);
X}
X
Xprint_bounds(box)
Xdouble box[2][3];
X{
X    extern FILE *fstats;
X    fprintf(fstats,"\tX: %f to %f\n",box[LOW][X], box[HIGH][X]);
X    fprintf(fstats,"\tY: %f to %f\n",box[LOW][Y], box[HIGH][Y]);
X    fprintf(fstats,"\tZ: %f to %f\n",box[LOW][Z], box[HIGH][Z]);
X}
END_OF_FILE
if test 5679 -ne `wc -c <'src/bounds.c'`; then
    echo shar: \"'src/bounds.c'\" unpacked with wrong size!
fi
# end of 'src/bounds.c'
fi
if test -f 'src/cone.c' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'src/cone.c'\"
else
echo shar: Extracting \"'src/cone.c'\" \(6421 characters\)
sed "s/^X//" >'src/cone.c' <<'END_OF_FILE'
X/*
X * cone.c
X *
X * Copyright (C) 1989, Craig E. Kolb
X *
X * This software may be freely copied, modified, and redistributed,
X * provided that this copyright notice is preserved on all copies.
X *
X * There is no warranty or other guarantee of fitness for this software,
X * it is provided solely .  Bug reports or fixes may be sent
X * to the author, who may or may not act on them as he desires.
X *
X * You may not include this software in a program or other software product
X * without supplying the source, or without informing the end-user that the
X * source is available for no extra charge.
X *
X * If you modify this software, you should include a notice giving the
X * name of the person performing the modification, the date of modification,
X * and the reason for such modification.
X *
X * $Id: cone.c,v 3.0 89/10/27 02:05:47 craig Exp $
X *
X * $Log:    cone.c,v $
X * Revision 3.0  89/10/27  02:05:47  craig
X * Baseline for first official release.
X *
X */
X#include <stdio.h>
X#include <math.h>
X#include "typedefs.h"
X#include "funcdefs.h"
X#include "constants.h"
X
XObject *
Xmakcone(surf, cent, ax, br, ar, trans)
Xchar *surf;
XVector *cent, *ax;
Xdouble br, ar;
XTransInfo *trans;
X{
X    Cone *cone;
X    Primitive *prim;
X    Object *newobj;
X    extern int yylineno, Quiet;
X    double len, dtmp;
X    Vector axis, base, tmp;
X
X    prim = mallocprim();
X    prim->surf = find_surface(surf);
X    prim->type = CONE;
X    newobj = new_object(NULL, CONE, (char *)prim, (Trans *)NULL);
X    cone = (Cone *)Malloc(sizeof(Cone));
X    prim->objpnt.p_cone = cone;
X
X    /*
X     * Cones are defined by a basepoint, an apex point, and
X     * base and apex radii.  The cone is stored as
X     * the origin of the cone, the change in radius per
X     * unit distance from the origin of the cone, the maximum z-value
X     * of the cone, and "start_pos",
X     * the distance along the axis from the cone origin where
X     * the first endcap appears (where the passed "basepoint"
X     * appears).
X     *
X     * The intcone() routine intersects a ray with a cone aligned
X     * along the Z axis.  Thus, we must define a transformation
X     * matrix which will transform an axis-aligned cone to the desired.
X     */
X
X    /*
X     * The passed basepoint must be closer to the origin of the
X     * cone than the apex point, implying that the base radius
X     * must be smaller than the apex radius.  If the values passed
X     * reflect the opposite, we switch everything.
X     */
X    if(ar < br) {
X        tmp = *cent;
X        *cent = *ax;
X        *ax = tmp;
X        dtmp = br;
X        br = ar;
X        ar = dtmp;
X    } else if (equal(ar, br)) {
X        /*
X         * If the base and apex radii are equal, then we
X         * can treat the cone as a cylinder.
X         */
X        return makcyl(surf, cent, ax, br, trans);
X    }
X    /*
X     * Find the axis and axis length.
X     */
X    vecsub(*ax, *cent, &axis);
X    len = normalize(&axis);
X    if (len < EPSILON) {
X        if (!Quiet)
X            fprintf(stderr,"Degenerate cone (line %d).\n",
X                            yylineno);
X        free((char *)cone);
X        free((char *)prim);
X        free((char *)newobj);
X        return (Object *)0;
X    }
X    cone->apex_rad = ar;
X    /*
X     * "tantheta" is the change in radius per unit length along
X     * the cone axis.
X     */
X    cone->tantheta = (ar - br) / len;
X    /*
X     * Start pos defines the point along the axis where the first
X     * endcap should be placed.
X     */
X    cone->start_pos = br / cone->tantheta;
X    /*
X     * Find the true base (origin) of the cone.
X     */
X    scalar_prod(-cone->start_pos, axis, &base);
X    vecadd(base, *cent, &base);
X    /*
X     * The apex endcap is placed cone->len units from the cone
X     * origin.
X     */
X    cone->end_pos = cone->start_pos + len;
X    /*
X     * Calculate rotation matrix to map from world space to cone space.
X     */
X/*    if (equal(axis.z*axis.z, 1.)) {
X        tmp.x = 0.;
X        tmp.y = -axis.z;
X        tmp.z = 0.;
X    } else { */
X        tmp.x = axis.y;
X        tmp.y = -axis.x;
X        tmp.z = 0.;
X    /*} */
X    rotate(trans, &tmp, acos(axis.z));
X    translate(trans, &base);
X    cone->tantheta *= cone->tantheta;
X
X    return newobj;
X}
X
X/*
X * Ray-cone intersection test.  This routine is far from optimal, but
X * it's straight-forward and it works...
X */
Xdouble
Xintcone(pos, ray, obj)
XVector           *pos, *ray;
XPrimitive       *obj;
X{
X    double t1, t2, a, b, c, disc, zpos, et1, et2;
X    double x, y;
X    extern unsigned long primtests[];
X    Cone *cone;
X
X    primtests[CONE]++;
X    cone = obj->objpnt.p_cone;
X
X    /*
X     * Recall that 'tantheta' is really tantheta^2.
X     */
X    a = ray->x * ray->x + ray->y * ray->y - ray->z*ray->z*cone->tantheta;
X    b = ray->x * pos->x + ray->y * pos->y - cone->tantheta*ray->z*pos->z;
X    c = pos->x*pos->x + pos->y*pos->y - cone->tantheta*pos->z*pos->z;
X
X    if (equal(a, 0.)) {
X        /*
X         * Only one intersection point...
X         */
X        t1 = -c / b;
X        zpos = pos->z + t1 * ray->z;
X        if (t1 < EPSILON || zpos < cone->start_pos ||
X            zpos > cone->end_pos)
X            t1 = FAR_AWAY;
X        t2 = FAR_AWAY;
X    } else {
X        disc = b*b - a*c;
X        if(disc < 0.)
X            return 0.;        /* No possible intersection */
X        disc = sqrt(disc);
X        t1 = (-b + disc) / a;
X        t2 = (-b - disc) / a;
X        /*
X         * Clip intersection points.
X         */
X        zpos = pos->z + t1 * ray->z;
X        if (t1 < EPSILON || zpos < cone->start_pos ||
X            zpos > cone->end_pos)
X            t1 = FAR_AWAY;
X        zpos = pos->z + t2 * ray->z;
X        if (t2 < EPSILON || zpos < cone->start_pos ||
X            zpos > cone->end_pos)
X            t2 = FAR_AWAY;
X    }
X    /*
X     * Find t for both endcaps.
X     */
X    et1 = (cone->start_pos - pos->z) / ray->z;
X    x = pos->x + et1 * ray->x;
X    y = pos->y + et1 * ray->y;
X    if (x*x + y*y > cone->start_pos*cone->start_pos*cone->tantheta)
X        et1 = FAR_AWAY;
X    et2 = (cone->end_pos - pos->z) / ray->z;
X    x = pos->x + et2 * ray->x;
X    y = pos->y + et2 * ray->y;
X    if (x*x + y*y > cone->end_pos*cone->end_pos*cone->tantheta)
X        et2 = FAR_AWAY;
X
X    t1 = min(t1, min(t2, min(et1, et2)));
X    return (t1 == FAR_AWAY ? 0. : t1);
X}
X
X/*
X * Compute the normal to a cone at a given location on its surface.
X */
Xnrmcone(pos, obj, nrm)
XVector *pos, *nrm;
XPrimitive *obj;
X{
X    Cone *cone;
X
X    cone = obj->objpnt.p_cone;
X
X    if (equal(pos->z, cone->start_pos)) {
X        nrm->x = nrm->y = 0.;
X        nrm->z = -1.;
X    } else if (equal(pos->z, cone->end_pos)) {
X        nrm->x = nrm->y = 0.;
X        nrm->z = 1.;
X    } else {
X        /*
X         * The following is equal to
X         * (pos X (0, 0, 1)) X pos
X         */
X        nrm->x = pos->x * pos->z;
X        nrm->y = pos->y * pos->z;
X        nrm->z = -pos->x * pos->x - pos->y * pos->y;
X    }
X}
X
X/*
X * Return the extent of a cone.
X */
Xconeextent(o, bounds)
XPrimitive *o;
Xdouble bounds[2][3];
X{
X    Cone *cone;
X
X    cone = o->objpnt.p_cone;
X
X    bounds[LOW][X] = bounds[LOW][Y] = -cone->apex_rad;
X    bounds[HIGH][X] = bounds[HIGH][Y] = cone->apex_rad;
X    bounds[LOW][Z] = cone->start_pos;
X    bounds[HIGH][Z] = cone->end_pos;
X}
END_OF_FILE
if test 6421 -ne `wc -c <'src/cone.c'`; then
    echo shar: \"'src/cone.c'\" unpacked with wrong size!
fi
# end of 'src/cone.c'
fi
if test -f 'src/cylinder.c' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'src/cylinder.c'\"
else
echo shar: Extracting \"'src/cylinder.c'\" \(4819 characters\)
sed "s/^X//" >'src/cylinder.c' <<'END_OF_FILE'
X/*
X * cylinder.c
X *
X * Copyright (C) 1989, Craig E. Kolb
X *
X * This software may be freely copied, modified, and redistributed,
X * provided that this copyright notice is preserved on all copies.
X *
X * There is no warranty or other guarantee of fitness for this software,
X * it is provided solely .  Bug reports or fixes may be sent
X * to the author, who may or may not act on them as he desires.
X *
X * You may not include this software in a program or other software product
X * without supplying the source, or without informing the end-user that the
X * source is available for no extra charge.
X *
X * If you modify this software, you should include a notice giving the
X * name of the person performing the modification, the date of modification,
X * and the reason for such modification.
X *
X * $Id: cylinder.c,v 3.0 89/10/27 02:05:48 craig Exp $
X *
X * $Log:    cylinder.c,v $
X * Revision 3.0  89/10/27  02:05:48  craig
X * Baseline for first official release.
X *
X */
X#include <stdio.h>
X#include <math.h>
X#include "typedefs.h"
X#include "funcdefs.h"
X#include "constants.h"
X
XObject *
Xmakcyl(surf, cent, ax, r, trans)
Xchar *surf;
XVector *cent, *ax;
Xdouble r;
XTransInfo *trans;
X{
X    Cylinder *cyl;
X    Primitive *prim;
X    Object *newobj;
X    double len;
X    extern int yylineno, Quiet;
X    Vector axis, dir;
X
X    if (r <= 0.) {
X        if (!Quiet)
X            fprintf(stderr,"Invalid cylinder radius (line %d)\n",
X                yylineno);
X        return (Object *)0;
X    }
X
X    prim = mallocprim();
X    newobj = new_object(NULL, CYL, (char *)prim, (Trans *)NULL);
X    prim->surf = find_surface(surf);
X    prim->type = CYL;
X    cyl = (Cylinder *)Malloc(sizeof(Cylinder));
X    prim->objpnt.p_cylinder = cyl;
X
X    axis.x = ax->x - cent->x;
X    axis.y = ax->y - cent->y;
X    axis.z = ax->z - cent->z;
X
X    len = normalize(&axis);
X    if(len < EPSILON) {
X        if (!Quiet)
X            fprintf(stderr,"Degenerate cylinder (line %d).\n",
X                            yylineno);
X        free((char *)cyl);
X        free((char *)prim);
X        return (Object *)0;
X    }
X
X    cyl->rad = r*r;
X    cyl->len = len;
X    /*
X     * Define matrix to transform from axis-aligned to desired cylinder.
X     */
X    dir.x = axis.y;
X    dir.y = -axis.x;
X    dir.z = 0.;
X    rotate(trans, &dir, acos(axis.z));
X    translate(trans, cent);
X
X    return newobj;
X}
X
X/*
X * Ray-cylinder intersection test.
X */
Xdouble
Xintcyl(pos, ray, obj)
XVector           *pos, *ray;
XPrimitive       *obj;
X{
X    double t1, t2, a, b, c, zpos1, zpos2, et1, et2, x, y, disc;
X    extern unsigned long primtests[];
X    Cylinder *cyl;
X
X    primtests[CYL]++;
X    cyl = obj->objpnt.p_cylinder;
X
X    a = ray->x * ray->x + ray->y * ray->y;
X    c = pos->x*pos->x + pos->y*pos->y - cyl->rad;
X
X    if (a < EPSILON*EPSILON) {    /* |ray->z| == 1. */
X        if(c < EPSILON*EPSILON)    /* Within endcap */
X            /* Wrong if origin is inside cylinder. */
X            return min(-pos->z / ray->z,
X                  (cyl->len - pos->z) / ray->z);
X        return 0.;
X    }
X
X    b = ray->x * pos->x + ray->y * pos->y;
X    disc = b*b - a*c;
X    if(disc < 0.)
X        return 0.;
X    disc = sqrt(disc);
X    t1 = (-b + disc) / a;
X    t2 = (-b - disc) / a;
X    if(t1 < EPSILON && t2 < EPSILON)
X        return 0.;
X    zpos1 = pos->z + t1 * ray->z;
X    zpos2 = pos->z + t2 * ray->z;
X    if ((zpos1 > cyl->len && zpos2 > cyl->len) ||
X        (zpos1 < 0. && zpos2 < 0.))
X        return 0.;
X    if (t1 < EPSILON)
X        t1 = FAR_AWAY;
X    if (t2 < EPSILON)
X        t2 = FAR_AWAY;
X    if (t1 == FAR_AWAY && t2 == FAR_AWAY)
X        return 0.;
X    /*
X     * Don't bother checking endcaps if both intersection points
X     * are on the cylinder.
X     */
X    if ((zpos1 > 0. && zpos1 < cyl->len && zpos2 > 0. && zpos2 < cyl->len))
X        return min(t1, t2);
X    /*
X      * It's possible to get rid of the ray-disc intersection tests
X     * (by looking at t1, t2 and zpos1, zpos), but the code gets messy.
X     */
X    if (zpos1 < 0. || zpos1 > cyl->len)
X        t1 = FAR_AWAY;
X    if (zpos2 < 0. || zpos2 > cyl->len)
X        t2 = FAR_AWAY;
X    et1 = -pos->z / ray->z;
X    x = pos->x + et1 * ray->x;
X    y = pos->y + et1 * ray->y;
X    if (x*x + y*y > cyl->rad)
X        et1 = FAR_AWAY;
X    et2 = (cyl->len - pos->z) / ray->z;
X    x = pos->x + et2 * ray->x;
X    y = pos->y + et2 * ray->y;
X    if (x*x + y*y > cyl->rad)
X        et2 = FAR_AWAY;
X    t1 = min(t1, min(t2, min(et1, et2)));
X    return (t1 == FAR_AWAY ? 0. : t1);
X}
X
Xnrmcyl(pos, obj, nrm)
XVector *pos, *nrm;
XPrimitive *obj;
X{
X    Cylinder *cyl;
X    double dist;
X
X    cyl = obj->objpnt.p_cylinder;
X
X    dist = pos->x*pos->x + pos->y*pos->y;
X    if (dist+EPSILON < cyl->rad) {
X        if (equal(pos->z,0.)) {
X            /*
X             * Hit on lower endcap.
X             */
X            nrm->x = nrm->y = 0.;
X            nrm->z = -1.;
X        } else {
X            /*
X             * Hit on upper endcap.
X             */
X            nrm->x = nrm->y = 0.;
X            nrm->z = 1.;
X        }
X    } else {    /* Hit along cylinder. */
X        nrm->x = pos->x;
X        nrm->y = pos->y;
X        nrm->z = 0.;
X        /* Will be normalized by ShadeRay(). */
X    }
X}
X
Xcylextent(o, bounds)
XPrimitive *o;
Xdouble bounds[2][3];
X{
X    Cylinder *cyl;
X    double r;
X    cyl = o->objpnt.p_cylinder;
X
X    r = sqrt(cyl->rad);
X    bounds[LOW][X] = bounds[LOW][Y] = -r;
X    bounds[HIGH][X] = bounds[HIGH][Y] = r;
X    bounds[LOW][Z] = 0.;
X    bounds[HIGH][Z] = cyl->len;
X}
END_OF_FILE
if test 4819 -ne `wc -c <'src/cylinder.c'`; then
    echo shar: \"'src/cylinder.c'\" unpacked with wrong size!
fi
# end of 'src/cylinder.c'
fi
if test -f 'src/intersect.c' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'src/intersect.c'\"
else
echo shar: Extracting \"'src/intersect.c'\" \(6514 characters\)
sed "s/^X//" >'src/intersect.c' <<'END_OF_FILE'
X/*
X * intersect.c
X *
X * Copyright (C) 1989, Craig E. Kolb
X *
X * This software may be freely copied, modified, and redistributed,
X * provided that this copyright notice is preserved on all copies.
X *
X * There is no warranty or other guarantee of fitness for this software,
X * it is provided solely .  Bug reports or fixes may be sent
X * to the author, who may or may not act on them as he desires.
X *
X * You may not include this software in a program or other software product
X * without supplying the source, or without informing the end-user that the
X * source is available for no extra charge.
X *
X * If you modify this software, you should include a notice giving the
X * name of the person performing the modification, the date of modification,
X * and the reason for such modification.
X *
X * $Id: intersect.c,v 3.0 89/10/27 02:05:53 craig Exp $
X *
X * $Log:    intersect.c,v $
X * Revision 3.0  89/10/27  02:05:53  craig
X * Baseline for first official release.
X *
X */
X#include <math.h>
X#include <stdio.h>
X#include "typedefs.h"
X#include "funcdefs.h"
X#include "constants.h"
X
X/*
X * Primitive intersection routines
X */
Xdouble (*objint[])() = {intsph, intbox, inttri, intsup, intplane, intcyl,
X            intpoly, inttri, intcone, inthf};
X/*
X * Primitive normal routines
X */
Xint (*objnrm[])() = {nrmsph, nrmbox, nrmtri, nrmsup, nrmplane, nrmcyl,
X            nrmpoly, nrmtri, nrmcone, nrmhf};
X/*
X * object extent box routines
X */
Xint (*objextent[])() = {sphextent, boxextent, triextent, supextent,
X            planeextent, cylextent, polyextent, triextent,
X            coneextent, hfextent};
X
Xunsigned long int primtests[PRIMTYPES], primhits[PRIMTYPES];
X
Xchar *primnames[PRIMTYPES] = {  "Sphere", "Box", "Triangle", "Superq", "Plane",
X                "Cylinder", "Polygon", "Phongtri", "Cone",
X                "Heightfield"};
X
X/*
X * Flags indicating whether or not we should check for intersection
X * with an object's bounding box before we check for intersection
X * with the object.
X */
Xchar CheckBounds[] = {TRUE, FALSE, TRUE, FALSE, FALSE, TRUE, TRUE, TRUE,
X            TRUE, FALSE, FALSE, FALSE};
X
X/*
X * Top-level raytracing routine.  Increment ray number, initialize
X * intersection information and trace ray through "World" object.
X */
Xdouble
XTraceRay(source, ray, hitinfo)
XPrimitive *source;
XRay *ray;
XHitInfo *hitinfo;
X{
X    extern Object *World;
X    extern double intersect();
X
X    return intersect(World, source, ray, hitinfo);
X}
X
X/*
X * Intersect object & ray.  Return distance from "pos" along "ray" to
X * intersection point.  Return value <= 0 indicates no intersection.
X */
Xdouble
Xintersect(obj, source, ray, hitinfo)
XObject *obj;                /* Object to be tested. */
XPrimitive *source;            /* Prim, if any, that pos is on. */
XRay *ray;                /* Ray origin, direction. */
XHitInfo *hitinfo;            /* Data on intersection (pos, norm) */
X{
X    Ray newray;
X    double dist, distfact, TransformRay();
X    extern int Cache;
X
X    /*
X     * Check ray/bounding volume intersection, if required.
X     */
X    if (CheckBounds[obj->type] &&
X        OutOfBounds(&ray->pos, obj->bounds) &&
X        IntBounds(ray, obj->bounds) < EPSILON)
X            return 0.;
X
X    newray = *ray;
X
X    /*
X     * Transform the ray if necessary.
X     */
X    if (obj->trans != (Trans *)0) {
X        /*
X         * Transforming the ray can change the distance between
X         * the ray origin and the point of intersection.
X         * We save the amount the ray is "stretched" and later
X         * divide the computed distance by this amount.
X         */
X        distfact = TransformRay(&newray, &obj->trans->world2obj);
X    }
X
X    /*
X     * Call correct intersection routine.
X     */
X    if (obj->type == GRID)
X        dist = int_grid((Grid *)obj->data, source, &newray, hitinfo);
X    else if (obj->type == LIST)
X        dist = int_list((List *)obj->data, source, &newray, hitinfo);
X    else
X        dist = int_primitive((Primitive *)obj->data, source, &newray,
X                    hitinfo);
X
X    if (dist < EPSILON)
X        return 0.;
X
X    /*
X     * If this is a shadow ray, don't bother with texture mapping
X     * or transformation of normal.
X     */
X    if (ray->shadow) {
X        if (obj->trans == (Trans *)0)
X            return dist;
X        else if (Cache)
X            /*
X              * Keep track of total transformation applied to ray
X              * if necessary.
X              */
X            mmult(hitinfo->totaltrans, &obj->trans->world2obj,
X                hitinfo->totaltrans);
X        return dist / distfact;
X    }
X    /*
X     * Perform texture mapping.
X     */
X    if (obj->texture)
X        apply_textures(hitinfo, obj->texture);
X
X    if (obj->trans) {
X        /*
X         * Transform hitinfo structure.  As things stand,
X         * this just means transforming the normal and
X         * dividing "dist" by the amount the ray was
X         * stretched.
X         */
X        dist /= distfact;
X        TransformNormal(&hitinfo->norm, &obj->trans->world2obj);
X    }
X
X    return dist;
X}
X
X/*
X * Intersect ray & primitive object.
X */
Xdouble
Xint_primitive(prim, source, ray, hitinfo)
XPrimitive *prim, *source;
XRay *ray;
XHitInfo *hitinfo;
X{
X    double dist;
X
X    if (prim == source && prim->type != HF)
X        /*
X         * Don't check for intersection with "source", unless
X         * source is a height field.  (Height fields may shadow
X         * themselves.)
X         */
X        return 0.;
X
X    dist = (*objint[prim->type]) (&ray->pos, &ray->dir, prim);
X
X    if (dist < EPSILON)
X        return 0.;
X
X    primhits[prim->type]++;
X    hitinfo->prim = prim;
X    hitinfo->surf = *prim->surf;
X
X    if (ray->shadow)
X        return dist;    /* If a shadow ray, don't bother with normal */
X    /*
X     * Calculate point of intersection in object space.
X      * (The point of intersection in world space is
X     * calculated in ShadeRay().)
X     */
X    addscaledvec(ray->pos, dist, ray->dir, &hitinfo->pos);
X
X    /*
X     * Find normal to primitive.
X     */
X    (*objnrm[prim->type]) (&hitinfo->pos, prim, &hitinfo->norm);
X
X    /*
X     * Make sure normal points towards ray origin.  If surface is
X     * transparent, keep as-is, as the normal indicates whether we're
X     * entering or exiting.  If the prim is a superquadric, don't flip,
X     * as this leads to strange edge effects.
X     */
X    if (dotp(&ray->dir, &hitinfo->norm) > 0 && hitinfo->surf.transp == 0. &&
X        prim->type != SUPERQ) {
X        scalar_prod(-1., hitinfo->norm, &hitinfo->norm);
X    }
X
X    return dist;
X}
X
Xprint_prim_stats()
X{
X    long int totaltests, totalhits;
X    extern FILE *fstats;
X    int i;
X
X    totaltests = totalhits = 0;
X    for (i = 0; i < PRIMTYPES; i++) {
X        if (primtests[i] == 0)
X            continue;
X        fprintf(fstats,"%s intersection tests:\t%ld (%ld hit, %f%%)\n",
X                primnames[i], primtests[i],
X                primhits[i],
X                100.*(float)primhits[i]/(float)primtests[i]);
X        totaltests += primtests[i];
X        totalhits += primhits[i];
X    }
X    fprintf(fstats,"Total intersection tests:\t%ld", totaltests);
X    if (totaltests == 0)
X        fprintf(fstats,"\n");
X    else
X        fprintf(fstats," (%ld hit, %f%%)\n", totalhits,
X            100.*(float)totalhits/(float)totaltests);
X}
END_OF_FILE
if test 6514 -ne `wc -c <'src/intersect.c'`; then
    echo shar: \"'src/intersect.c'\" unpacked with wrong size!
fi
# end of 'src/intersect.c'
fi
if test -f 'src/light.c' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'src/light.c'\"
else
echo shar: Extracting \"'src/light.c'\" \(5124 characters\)
sed "s/^X//" >'src/light.c' <<'END_OF_FILE'
X/*
X * light.c
X *
X * Copyright (C) 1989, Craig E. Kolb
X *
X * This software may be freely copied, modified, and redistributed,
X * provided that this copyright notice is preserved on all copies.
X *
X * There is no warranty or other guarantee of fitness for this software,
X * it is provided solely .  Bug reports or fixes may be sent
X * to the author, who may or may not act on them as he desires.
X *
X * You may not include this software in a program or other software product
X * without supplying the source, or without informing the end-user that the
X * source is available for no extra charge.
X *
X * If you modify this software, you should include a notice giving the
X * name of the person performing the modification, the date of modification,
X * and the reason for such modification.
X *
X * $Id: light.c,v 3.0 89/10/27 16:17:22 craig Exp $
X *
X * $Log:    light.c,v $
X * Revision 3.0  89/10/27  16:17:22  craig
X * Baseline for first official release.
X *
X */
X#include <stdio.h>
X#include <math.h>
X#include "typedefs.h"
X#include "funcdefs.h"
X#include "constants.h"
X
Xint    nlight;            /* # of lights defined */
Xint    NoShadows;        /* Don't trace shadow rays */
Xint    Cache = TRUE;        /* Use shadow-caching */
Xint    ClearShadows;        /* Shadow rays pass through transp. objects */
XLight    light[LIGHTS];        /* array of lights */
Xdouble    lightdist;        /* distance to light */
Xunsigned long CacheWorked, CacheFailed, ShadowHits;
X/*
X * Calculate ray from position to light # lnum.
X */
Xlightray(lp, objpos, lray)
XLight *lp;
XVector *objpos, *lray;
X{
X    if(lp->type == DIRECTIONAL) {
X        /*
X         * Directional sources only have direction.
X         */
X        *lray = lp->pos;
X        lightdist = FAR_AWAY;
X    } else {
X        /*
X         * Calculate ray from position to center of
X         * light source.
X         */
X        vecsub(lp->pos, *objpos, lray);
X        lightdist = normalize(lray);
X    }
X}
X
X/*
X * Find a coordinate system perpendicular to the ray from
X * a point of intersection to the center of an extended light source.
X */
XLightCoordSys(lp, pos, vector, xaxis, yaxis)
XLight *lp;
XVector *pos, *vector, *xaxis, *yaxis;
X{
X    /*
X     * Vector should *not* be normalized, xaxis & yaxis should be.
X     */
X    vector->x = lp->pos.x - pos->x;
X    vector->y = lp->pos.y - pos->y;
X    vector->z = lp->pos.z - pos->z;
X    xaxis->x = vector->y;
X    xaxis->y = -vector->x;
X    xaxis->z = 0.;
X    if (normalize(xaxis) == 0.) {
X        xaxis->x = 0.;
X        xaxis->y = -vector->z;
X        xaxis->z = vector->y;
X        if (normalize(xaxis) == 0.)
X            fprintf(stderr,"LightCoordSys: Can't find X axis!\n");
X        yaxis->x = (vector->y * vector->y) + (vector->z * vector->z);
X        yaxis->y = -vector->x * vector->y;
X        yaxis->z = -vector->x * vector->z;
X        (void)normalize(yaxis);
X    } else {
X        yaxis->x = vector->x * vector->z;
X        yaxis->y = vector->y * vector->z;
X        yaxis->z = -(vector->x * vector->x) -(vector->y * vector->y);
X        (void)normalize(yaxis);
X    }
X}
X
X/*
X * Trace ray from point of intersection to a light.  If an intersection
X * occurs at a distance less than "lightdist" (the distance to the
X * light source), then the point is in shadow, and 0 is returned.
X * Otherwise, the brightness (color) of the light is returned.  This
X * color may be modulated by any translucent objects which fall between
X * the point of intersection and the light source.
X */
Xinshadow(result, source, lp, pos, ray)
XColor *result;
XPrimitive *source;
XLight *lp;
XVector *pos, *ray;
X{
X    double s;
X    Ray tmpray, tray;
X    HitInfo hitinfo;
X    double atten, totaldist, TransformRay();
X    extern int level;
X    extern unsigned long ShadowRays;
X    extern double TraceRay();
X
X    if (NoShadows) {
X        *result = lp->color;
X        return FALSE;
X    }
X
X    ShadowRays++;
X    tmpray.pos = *pos;
X    tmpray.dir = *ray;    /* Medium not needed. */
X    tmpray.shadow = TRUE;
X    hitinfo.totaltrans = &lp->trans[level];
X
X    /*
X     * Check shadow cache if necessary.  (The following implies
X     * ... && Cache)
X     */
X    if (lp->cache[level]) {
X        tray = tmpray;
X        s = TransformRay(&tray, &lp->trans[level]);
X        s = int_primitive(lp->cache[level], source, &tray, &hitinfo)/s;
X        if (s > EPSILON && s < lightdist) {
X            CacheWorked++;
X            return TRUE;
X        }
X        CacheFailed++;
X        lp->cache[level] = (Primitive *)0;
X    }
X
X    if (Cache)
X        init_trans(hitinfo.totaltrans);
X
X    s = TraceRay(source, &tmpray, &hitinfo);
X
X    if (s < EPSILON || s > lightdist) {
X        *result = lp->color;
X        return FALSE;        /* Not in shadow. */
X    }
X
X    /*
X     * Otherwise, we've hit something.
X     */
X    ShadowHits++;
X    if (!ClearShadows || hitinfo.surf.transp == 0.) {
X        if (Cache)
X            lp->cache[level] = hitinfo.prim;
X        return TRUE;
X    }
X    /*
X     * We've hit a transparent object.  Attenuate the color
X     * of the light source and continue the ray until
X     * we hit background or a non-transparent object.
X     * Note that this is incorrect if any of the surfaces hit (or
X     * DefIndex) have differing indices of refraction.
X     */
X    atten = 1.;
X    totaldist = s;
X    do {
X        atten *= hitinfo.surf.transp;
X        if (atten < EPSILON)
X            return TRUE;
X        addscaledvec(tmpray.pos, s, tmpray.dir, &tmpray.pos);
X        /*
X         * Trace ray starting at new origin and in the
X         * same direction.
X         */
X        s = TraceRay(hitinfo.prim, &tmpray, &hitinfo);
X        totaldist += s;
X    } while (s > EPSILON && totaldist < lightdist);
X
X    ScaleColor(atten, lp->color, result);
X    return FALSE;
X}
END_OF_FILE
if test 5124 -ne `wc -c <'src/light.c'`; then
    echo shar: \"'src/light.c'\" unpacked with wrong size!
fi
# end of 'src/light.c'
fi
if test -f 'src/object.c' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'src/object.c'\"
else
echo shar: Extracting \"'src/object.c'\" \(4915 characters\)
sed "s/^X//" >'src/object.c' <<'END_OF_FILE'
X/*
X * object.c
X *
X * Copyright (C) 1989, Craig E. Kolb
X *
X * This software may be freely copied, modified, and redistributed,
X * provided that this copyright notice is preserved on all copies.
X *
X * There is no warranty or other guarantee of fitness for this software,
X * it is provided solely .  Bug reports or fixes may be sent
X * to the author, who may or may not act on them as he desires.
X *
X * You may not include this software in a program or other software product
X * without supplying the source, or without informing the end-user that the
X * source is available for no extra charge.
X *
X * If you modify this software, you should include a notice giving the
X * name of the person performing the modification, the date of modification,
X * and the reason for such modification.
X *
X * $Id: object.c,v 3.0 89/10/27 02:05:58 craig Exp $
X *
X * $Log:    object.c,v $
X * Revision 3.0  89/10/27  02:05:58  craig
X * Baseline for first official release.
X *
X */
X#include <stdio.h>
X#include <math.h>
X#include "constants.h"
X#include "typedefs.h"
X#include "funcdefs.h"
X#include "texture.h"
X
XObject *World;            /* World Object */
XObjList *Objects;        /* Linked list of defined objects */
Xint WorldXSize, WorldYSize, WorldZSize;        /* World grid resolution */
X
X/*
X * Create a new object with the given properties.
X */
XObject *
Xnew_object(name, type, data, trans)
Xchar *name, *data;
Xchar type;
XTrans *trans;
X{
X    Object *new;
X
X    new = (Object *)share_malloc(sizeof(Object));
X    new->name = strsave(name);
X    new->type = type;
X    new->data = data;
X    new->trans = trans;
X#ifdef LINDA
X    /*
X     * If the counter is in shared memory, processes will
X     * be modifying it left-and-right.  So, we cheat and
X     * make counter a pointer to a non-shared location and
X     * store the value there.
X     */
X    new->counter = (unsigned long *)malloc(sizeof(unsigned long));
X    *new->counter = 0;
X#else
X    new->counter = 0;
X#endif
X    new->texture = (Texture *)0;
X    /*
X     * bounds is left uninitialized.
X     */
X    return new;
X}
X
X/*
X * Add a copy of the named object to parent object
X */
XObject *
Xadd_child_named(name, parent)
Xchar *name;
XObject *parent;
X{
X    Object *child, *newobj;
X    int i;
X
X    child = get_object_named(name);
X    if (child == (Object *)0) {
X        fprintf(stderr,"There is no object named \"%s\".\n",name);
X        exit(1);
X    }
X    /*
X     * Create new object that points to child
X     * and add to 'parent' list.
X     */
X    newobj = add_child(child, parent);
X    /*
X     * New object's bounding box is initally the same
X     * as the child's.
X     */
X    for (i = 0; i < 3; i++) {
X        newobj->bounds[0][i] = child->bounds[0][i];
X        newobj->bounds[1][i] = child->bounds[1][i];
X    }
X    return newobj;
X}
X
X/*
X * Add primitive object to parent object.
X */
Xadd_prim(child, parent)
XObject *child, *parent;
X{
X    ObjList *newnode;
X
X    newnode = (ObjList *)share_malloc(sizeof(ObjList));
X    newnode->data = child;
X    if (parent == (Object *)0) {
X        newnode->next = (ObjList *)World->data;
X        World->data = (char *)newnode;
X    } else {
X        newnode->next = (ObjList *)parent->data;
X        parent->data = (char *)newnode;
X    }
X}
X
X/*
X * Make a copy of "child" and attach it to parent's linked list
X * of objects.
X */
XObject *
Xadd_child(child, parent)
XObject *child, *parent;
X{
X    Object *newobj;
X    ObjList *newnode;
X
X    newobj = new_object(NULL, child->type, child->data, child->trans);
X    newobj->texture = child->texture;
X    newnode = (ObjList *)share_malloc(sizeof(ObjList));
X    newnode->data = newobj;
X    if (parent == (Object *)0) {
X        newnode->next = (ObjList *)World->data;
X        World->data = (char *)newnode;
X    } else {
X        newnode->next = (ObjList *)parent->data;
X        parent->data = (char *)newnode;
X    }
X    return newobj;
X}
X
X/*
X * Return pointer to named object, NULL if no such object has been defined.
X */
XObject *
Xget_object_named(name)
Xchar *name;
X{
X    ObjList *ltmp;
X    for (ltmp = Objects; ltmp; ltmp = ltmp->next)
X        if (strcmp(name, ltmp->data->name) == 0)
X            return ltmp->data;
X    return (Object *)0;
X}
X
X/*
X * Add object to list of defined objects.
X */
Xadd_to_objects(obj)
XObject *obj;
X{
X    ObjList *ltmp;
X    extern int Verbose;
X    extern FILE *fstats;
X
X    ltmp = (ObjList *)Malloc(sizeof(ObjList));
X    ltmp->data = obj;
X    ltmp->next = Objects;
X    Objects = ltmp;
X    if (Verbose) {
X        /*
X         * Report bounding box of named object.
X         */
X        fprintf(fstats,"Object \"%s\" extent:\n", obj->name);
X        print_bounds(obj->bounds);
X    }
X}
X
X/*
X * Allocate space for a string, copy string into space.
X */
Xchar *
Xstrsave(s)
Xchar *s;
X{
X    char *tmp;
X
X    if (s == (char *)0)
X        return (char *)0;
X
X    tmp = (char *)Malloc((unsigned)strlen(s) + 1);
X    strcpy(tmp, s);
X    return tmp;
X}
X
X/*
X * Set "bounds" of primitive to be the extent of the primitive.
X */
Xset_prim_bounds(obj)
XObject *obj;
X{
X    extern int (*objextent[])();
X
X    (*objextent[((Primitive *)obj->data)->type])
X        ((Primitive *)obj->data, obj->bounds);
X    obj->bounds[LOW][X] -= EPSILON;
X    obj->bounds[HIGH][X] += EPSILON;
X    obj->bounds[LOW][Y] -= EPSILON;
X    obj->bounds[HIGH][Y] += EPSILON;
X    obj->bounds[LOW][Z] -= EPSILON;
X    obj->bounds[HIGH][Z] += EPSILON;
X}
X
END_OF_FILE
if test 4915 -ne `wc -c <'src/object.c'`; then
    echo shar: \"'src/object.c'\" unpacked with wrong size!
fi
# end of 'src/object.c'
fi
if test -f 'src/ray_options.c' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'src/ray_options.c'\"
else
echo shar: Extracting \"'src/ray_options.c'\" \(6048 characters\)
sed "s/^X//" >'src/ray_options.c' <<'END_OF_FILE'
X/*
X * ray_options.c
X *
X * Copyright (C) 1989, Craig E. Kolb
X *
X * This software may be freely copied, modified, and redistributed,
X * provided that this copyright notice is preserved on all copies.
X *
X * There is no warranty or other guarantee of fitness for this software,
X * it is provided solely .  Bug reports or fixes may be sent
X * to the author, who may or may not act on them as he desires.
X *
X * You may not include this software in a program or other software product
X * without supplying the source, or without informing the end-user that the
X * source is available for no extra charge.
X *
X * If you modify this software, you should include a notice giving the
X * name of the person performing the modification, the date of modification,
X * and the reason for such modification.
X *
X * $Id: ray_options.c,v 3.0 89/10/27 02:06:00 craig Exp $
X *
X * $Log:    ray_options.c,v $
X * Revision 3.0  89/10/27  02:06:00  craig
X * Baseline for first official release.
X *
X */
X#include <stdio.h>
X#ifdef SYSV
X#include <string.h>
X#else
X#include <strings.h>
X#endif
X#include "constants.h"
X#include "typedefs.h"
X
Xint Verbose;            /* Blabbering flag */
Xint TrashBadPoly;        /* Discard even mildly bad polygons */
Xint Quiet;            /* Don't be so verbose flag */
Xchar *progname;            /* argv[0] */
X
X#ifdef LINDA
Xint Workers;            /* # of workers */
X#endif
X
Xparse_options(argc, argv)
Xint argc;
Xchar **argv;
X{
X    extern char *infilename;
X    extern double RedContrast, GreenContrast, BlueContrast, atof();
X    extern int pixel_div, JitSamples, Xres, Yres, Jittered, Cache;
X    extern int Stereo, StartLine, Appending, NoShadows, ClearShadows;
X    extern double Separation, TreeCutoff;
X    extern char outfilename[];
X    extern FILE *fstats;
X
X    progname = argv[0];
X    fstats = stderr;
X
X    while(--argc) {
X        argv++;
X        if(argv[0][0] != '-')
X            break;
X        switch(argv[0][1]) {
X            case 'C':
X                RedContrast = atof(argv[1]);
X                GreenContrast = atof(argv[2]);
X                BlueContrast = atof(argv[3]);
X                argv += 3;
X                argc -= 3;
X                break;
X            case 'c':
X                ClearShadows = TRUE;
X                break;
X            case 'E':
X                Separation = atof(argv[1]);
X                argc--; argv++;
X                break;
X            case 'h':
X                usage();
X                exit(0);
X                break;
X            case 'j':
X                Jittered = TRUE;
X                break;
X            case 'L':
X                StartLine = atoi(argv[1]);
X                Appending = TRUE;
X                argc--; argv++;
X                break;
X            case 'l':
X                Stereo = LEFT;
X                break;
X            case 'n':
X                NoShadows = TRUE;
X                break;
X            case 'O':
X                strcpy(outfilename, argv[1]);
X                argv++;
X                argc--;
X                break;
X            case 'P':
X                pixel_div = atoi(argv[1]);
X                if(pixel_div < 0)
X                    pixel_div = 0;
X                argv++;
X                argc--;
X                break;
X            case 'p':
X                TrashBadPoly = TRUE;
X                break;
X            case 'q':
X                Quiet = TRUE;
X                break;
X            case 'R':
X                Xres = atoi(argv[1]);
X                Yres = atoi(argv[2]);
X                argv += 2;
X                argc -= 2;
X                break;
X            case 'r':
X                Stereo = RIGHT;
X                break;
X            case 'S':
X                JitSamples = atoi(argv[1]);
X                if (JitSamples < 1)
X                    JitSamples = 1;
X                argv++; argc--;
X                break;
X            case 's':
X                Cache = !Cache;
X                break;
X            case 'T':
X                TreeCutoff = atof(argv[1]);
X                argv++; argc--;
X                break;
X            case 'v':
X                Verbose = TRUE;
X                break;
X            case 'V':
X                Verbose = TRUE;
X                if (argv[1][0] == '-') {
X                    /* User probably blew it, and
X                     * it's difficult to remove a file
X                     * that begins with '-'...
X                     */
X                    usage();
X                    exit(2);
X                }
X                fstats = fopen(argv[1], "w");
X                if (fstats == (FILE *)0) {
X                    fprintf(stderr,"Cannot write to stats file %s\n",argv[0]);
X                    exit(2);
X                }
X                argv++; argc--;
X                break;
X#ifdef LINDA
X            case 'W':
X                Workers = atoi(argv[1]);
X                if (Workers < 0 || Workers > 17) {
X                    fprintf(stderr,"%d workers?!?\n",
X                            Workers);
X                    exit(3);
X                }
X                argv++; argc--;
X                break;
X#endif
X            default:
X                fprintf(stderr,"Bad argument: \"%s\"\n",argv[0]);
X                usage();
X                exit(1);
X        }
X    }
X
X    if(argc > 1) {
X        usage();
X        exit(1);
X    } else if(argc == 1)
X        infilename = argv[0];
X    else
X        infilename = (char *)NULL;
X
X    /*
X     * Although the user may have defined the output file name
X     * in the input file, it's best to force them to give a filename
X     * on the command line if using the -L option.  This saves situations
X     * where they forget to specify an output file to append to but
X     * aren't informed of it until startpic() is called (after the
X     * entire input file is read).
X     */
X    if (Appending && *outfilename == (char)NULL) {
X        fprintf(stderr,"The -L option requires the -O option.\n");
X        exit(4);
X    }
X
X    if (Stereo && Separation == 0.) {
X        fprintf(stderr,"You must specify eye separation (-E) in order ");
X        fprintf(stderr,"to enable Stereo mode.\n");
X        exit(4);
X    }
X}
X
Xusage()
X{
X    fprintf(stderr,"usage: %s [options] [filename]\n", progname);
X    fprintf(stderr,"Where options include:\n");
X    fprintf(stderr,"\t-C r g b\t(Set contrast threshold (0. - 1.).)\n");
X    fprintf(stderr,"\t-c \t\t(Trace shadow rays through clear objects.)\n");
X    fprintf(stderr,"\t-E eye_sep\t(Set eye separation.)\n");
X    fprintf(stderr,"\t-h \t\t(Print this message.)\n");
X    fprintf(stderr,"\t-j \t\t(Antialias using jittered sampling.)\n");
X    fprintf(stderr,"\t-L line#\t(Begin rendering at specified line.)\n");
X    fprintf(stderr,"\t-l \t\t(Render image for left eye view.)\n");
X    fprintf(stderr,"\t-n \t\t(Don't compute shadows.)\n");
X    fprintf(stderr,"\t-O outfile \t(Specify output file name.)\n");
X    fprintf(stderr,"\t-P pixel_divs\t(Set max depth for adaptive supersampling.)\n");
X    fprintf(stderr,"\t-p \t\t(Discard polygons with degenerate edges.)\n");
X    fprintf(stderr,"\t-q \t\t(Run quietly.)\n");
X    fprintf(stderr,"\t-R xres yres\t(Render at given resolution.)\n");
X    fprintf(stderr,"\t-r \t\t(Render image for right eye view.)\n");
X    fprintf(stderr,"\t-S samples\t(Use samples^2 jittered samples.)\n");
X    fprintf(stderr,"\t-s \t\t(Don't cache shadowing information.)\n");
X    fprintf(stderr,"\t-T thresh\t(Set adaptive ray tree cutoff value.)\n");
X    fprintf(stderr,"\t-V filename \t(Write verbose output to filename.)\n");
X    fprintf(stderr,"\t-v \t\t(Verbose output.)\n");
X#ifdef LINDA
X    fprintf(stderr,"\t-W workers (Specify number of worker processes.)\n");
X#endif
X}
END_OF_FILE
if test 6048 -ne `wc -c <'src/ray_options.c'`; then
    echo shar: \"'src/ray_options.c'\" unpacked with wrong size!
fi
# end of 'src/ray_options.c'
fi
if test -f 'src/voxels.c' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'src/voxels.c'\"
else
echo shar: Extracting \"'src/voxels.c'\" \(4916 characters\)
sed "s/^X//" >'src/voxels.c' <<'END_OF_FILE'
X/*
X * voxels.c
X *
X * Copyright (C) 1989, Craig E. Kolb
X *
X * This software may be freely copied, modified, and redistributed,
X * provided that this copyright notice is preserved on all copies.
X *
X * There is no warranty or other guarantee of fitness for this software,
X * it is provided solely .  Bug reports or fixes may be sent
X * to the author, who may or may not act on them as he desires.
X *
X * You may not include this software in a program or other software product
X * without supplying the source, or without informing the end-user that the
X * source is available for no extra charge.
X *
X * If you modify this software, you should include a notice giving the
X * name of the person performing the modification, the date of modification,
X * and the reason for such modification.
X *
X * $Id: voxels.c,v 3.0 89/10/27 02:06:09 craig Exp $
X *
X * $Log:    voxels.c,v $
X * Revision 3.0  89/10/27  02:06:09  craig
X * Baseline for first official release.
X *
X */
X#include <math.h>
X#include <stdio.h>
X#include "constants.h"
X#include "typedefs.h"
X#include "funcdefs.h"
X
X/*
X * Process World object, converting to a Grid or List.
X */
XSetupWorld()
X{
X    extern FILE *fstats;
X    extern Object *World;
X    extern int WorldXSize, WorldYSize, WorldZSize, Verbose;
X
X    if (World->type == GRID)
X        list2grid(World, WorldXSize, WorldYSize, WorldZSize);
X    else
X        make_list(World);
X
X    if (Verbose) {
X        fprintf(fstats,"World extent:\n");
X        print_bounds(World->bounds);
X    }
X}
X
X/*
X * Add object to grid's unbounded list.
X */
Xmake_unbounded(obj, grid)
XObject *obj;
XGrid *grid;
X{
X    ObjList *tmp;
X
X    tmp = (ObjList *)Malloc(sizeof(ObjList));
X
X    tmp->data = obj;
X    tmp->next = grid->unbounded;
X    grid->unbounded = tmp;
X}
X
X/*
X * Place an object in a grid.
X */
Xengrid(obj, grid)
XObject *obj;
XGrid *grid;
X{
X    int x, y, z, low[3], high[3];
X    ObjList *ltmp;
X
X    /*
X     * This routine should *never* be passed an unbounded object, but...
X     */
X    if (pos2grid(grid, obj->bounds[LOW], low) == 0 ||
X        pos2grid(grid, obj->bounds[HIGH], high) == 0 ||
X        obj->bounds[LOW][X] > obj->bounds[HIGH][X]) {
X        /*
X         * Object is partially on wholly outside of
X         * grid -- this should never happen, but just
X         * in case...
X         */
X        make_unbounded(obj, grid);
X        fprintf(stderr,"Strange, engrid got an unbounded object...\n");
X        return;
X        }
X
X    /*
X     * For each voxel that intersects the object's bounding
X     * box, add pointer to this object to voxel's linked list.
X      */
X    for (x = low[X]; x <= high[X]; x++) {
X        for (y = low[Y]; y <= high[Y]; y++) {
X            for (z = low[Z]; z <= high[Z]; z++) {
X                ltmp = (ObjList *)share_malloc(sizeof(ObjList));
X                ltmp->data = obj;
X                ltmp->next = grid->cells[x][y][z];
X                grid->cells[x][y][z] = ltmp;
X            }
X        }
X    }
X}
X
X/*
X * Convert 3D point to index into grid's voxels.
X */
Xpos2grid(grid, pos, index)
XGrid *grid;
Xdouble pos[3];
Xint index[3];
X{
X    index[X] = (int)(x2voxel(grid, pos[0]));
X    index[Y] = (int)(y2voxel(grid, pos[1]));
X    index[Z] = (int)(z2voxel(grid, pos[2]));
X
X    if (index[X] == grid->xsize)
X        index[X] = grid->xsize -1;
X    if (index[Y] == grid->ysize)
X        index[Y] = grid->ysize -1;
X    if (index[Z] == grid->zsize)
X        index[Z] = grid->zsize -1;
X
X    if (index[X] < 0 || index[X] >= grid->xsize ||
X        index[Y] < 0 || index[Y] >= grid->ysize ||
X        index[Z] < 0 || index[Z] >= grid->zsize)
X        return 0;
X    return 1;
X}
X
X/*
X * Convert a linked list of objects to a Grid.
X */
Xlist2grid(obj, xsize, ysize, zsize)
XObject *obj;
Xint xsize, ysize, zsize;
X{
X    Grid *grid;
X    Object *otmp;
X    ObjList *ltmp;
X    int x, y, i;
X    extern ObjList *find_bounds();
X
X    grid = (Grid *)Malloc(sizeof(Grid));
X
X    /*
X     * Find bounding box of bounded objects and get list of
X     * unbounded objects.
X     */
X    grid->unbounded = find_bounds((ObjList **)&obj->data, obj->bounds);
X
X    grid->xsize = xsize; grid->ysize = ysize; grid->zsize = zsize;
X
X    for (i = 0; i < 3; i++) {
X        obj->bounds[LOW][i] -= 2. * EPSILON;
X        obj->bounds[HIGH][i] += 2. * EPSILON;
X        grid->bounds[LOW][i] = obj->bounds[LOW][i];
X        grid->bounds[HIGH][i] = obj->bounds[HIGH][i];
X    }
X    grid->voxsize[X] = (grid->bounds[HIGH][X]-grid->bounds[LOW][X])/xsize;
X    grid->voxsize[Y] = (grid->bounds[HIGH][Y]-grid->bounds[LOW][Y])/ysize;
X    grid->voxsize[Z] = (grid->bounds[HIGH][Z]-grid->bounds[LOW][Z])/zsize;
X
X    /*
X     * Allocate voxels.
X     */
X    grid->cells = (ObjList ****)share_malloc(xsize * sizeof(ObjList ***));
X    for (x = 0; x < xsize; x++) {
X        grid->cells[x] = (ObjList ***)share_malloc(ysize*sizeof(ObjList **));
X        for (y = 0; y < ysize; y++)
X            grid->cells[x][y] = (ObjList **)share_calloc((unsigned)zsize,
X                            sizeof(ObjList *));
X    }
X
X    /*
X     * obj->data now holds a linked list of bounded objects.
X     */
X    for(ltmp = (ObjList *)obj->data; ltmp; ltmp = ltmp->next) {
X        otmp = ltmp->data;
X        engrid(otmp, grid);
X        free(ltmp);
X    }
X    obj->type = GRID;
X    obj->data = (char *)grid;
X}
X#ifdef MULTIMAX
X
Xchar *
Xshare_calloc(num, siz)
Xint num;
Xunsigned int siz;
X{
X    char *res;
X
X    res = share_malloc(num*siz);
X    bzero(res, num*siz);
X    return res;
X}
X#endif
END_OF_FILE
if test 4916 -ne `wc -c <'src/voxels.c'`; then
    echo shar: \"'src/voxels.c'\" unpacked with wrong size!
fi
# end of 'src/voxels.c'
fi
echo shar: End of archive 3 \(of 8\).
cp /dev/null ark3isdone
MISSING=""
for I in 1 2 3 4 5 6 7 8 ; do
    if test ! -f ark${I}isdone ; then
    MISSING="${MISSING} ${I}"
    fi
done
if test "${MISSING}" = "" ; then
    echo You have unpacked all 8 archives.
    rm -f ark[1-9]isdone
else
    echo You still need to unpack the following archives:
    echo "        " ${MISSING}
fi
##  End of shell archive.
exit 0